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Abstract 

Along the lines of the nonlinear response theory developed by Ruelle, in a previous paper we 
have proved under rather general conditions that Kramers-Kronig dispersion relations and sum 
rules apply for a class of susceptibilities describing at any order of perturbation the response of 
Axiom A non equilibrium steady state systems to weak monochromatic forcings. We present here 
the first evidence of the validity of these integral relations for the linear and the second harmonic 
response for the perturbed Lorenz 63 system, by showing that numerical simulations agree up to 
high degree of accuracy with the theoretical predictions. Some new theoretical results, showing 
how to obtain recursively harmonic generation susceptibilities for general observables, are also 
presented. Our findings confirm the conceptual validity of the nonlinear response theory, suggest 
that the theory can be extended for more general non equilibrium steady state systems, and shed 
new light on the applicability of very general tools, based only upon the principle of causality, for 
diagnosing the behavior of perturbed chaotic systems and reconstructing their output signals, in 
situations where the fluctuation-dissipation relation is not of great help. 
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I. INTRODUCTION 



A rather wide class of physical problems can be framed as the analysis of the sensitivity 
of the statistical properties of a system to external parameters, and, actually, such sensitiv- 
ities often define physical quantities of great conceptual relevance, as in outstanding case of 
Maxwell's relations in classical thermodynamics. In this context, the response theory for- 
malizes the thought experimental apparatus comprising of the system under investigation, 
of a measuring devic, of a clock, and of a set of turnable knobs controlling, possibly with 
continuity, the value of the external parameters. 

In the case of physical systems near to an equilibrium represented by the canonical mea- 
sure, Kubo introduced a response theory able to describe up to any order of nonlinearity 
the impact of weak perturbations to the Hamiltonian on the statistical properties of a gen- 
eral observable. Apart from its profound theoretical value, the Kubo approach has allowed 
for the definition of applicable tools for obtaining a systematic and detailed analysis of the 
forced fluctuations due to weak external forcings, thus proving to be a suitable tool for a 
comprehensive treatment of optical, acoustical, mechanical phenomena {4]. 

The response theory recently developed by Ruelle [l, 2] for the description of the impact 
of small, periodic perturbations on the statistical properties of a measurable observable of 
an Axiom A dynamical system - whose attractor defines an SRB measure - allows for 
a generalization of the whole apparatus of the Kubo theory to the case of non-equilibrium 
steady state statistical mechanical systems. The Kubo and Ruelle response theories are 
formally rather similar, as in both cases the impact of perturbation can be expressed as 
the expectation value of a specific operator over the statistical ensemble of the unperturbed 
state, the main (crucial) diversity being, instead, in the properties of the probability measure 
of integration. 

In quasi-equilibrium statistical mechanics, it has been recognized that, thanks to the gen- 
eral principle of causality of the response, the linear and nonlinear generalized susceptibilities 
of the system, describing its response in the frequency domain, feature specific properties of 
analyticity, allowing for the definition of Kramers-Kronig (K-K) relations, connecting their 
real and imaginary parts, and for the derivation of sum rules, which connect the response of 
the system at all 'requencies to the expectation value of specific observables at the unper- 



turbed state [Sl, lid, lid, llSl, l22| . These integral relations play a crucial role in several research 



and technological areas, and most prominently in condensed matter physics and material 
science [la, ■ 



Ruelle [l| proved that, whereas the fluctuation-dissipation relation, cornerstone of 
quasi-equilibrium statistical mechanics, cannot be straightforwardly extended to the non- 
equilibrium case, thus clarifying an earlier intuition by Lorenz on the non equivalence be- 
tween forced a„d .ee flu.ua.on. i,UKe curate .yste. H, the case of perturbed Axiom 
A system it is possible to define a linear susceptibility which is analytic in the upper complex 
plane of the frequency variable, and derive rigorously K-K relations. More recently, building 



upon these results and following [l^, in 2j] it was shown that a large and well de- 
fined class of nonlinear susceptibilities, including those responsible for harmonic generation 
processes at all orders, obey generalized K-K relations and, additionally, sum rules can be 
established. 

It is reasonable to expect that a much larger set of dynamical systems than the Axiom 
A ones obeys these properties, as recently suggested by theoretical argunrents Q and by 
some rather convincing numerical simulations [17|, 119|. Note that, nevertheless, results on 
Axiom A systems are already quite powerful in terms of practical applications, if one accepts 
the chaotic hypothesis by Gallavotti and Cohen [sl, stating that large systems behave as 
though they were Axiom A systems when macroscopic statistical properties are considered. 

In this paper we present the first evidence of applicability of dispersion relations and 
sum rules for a nonlinear susceptibility of a chaotic dynamical system. We consider a 
monochromatic (weak)j external perturbation to the flow of the Lorenz 63 system [l^. First, 
we perform the analysis of the linear response, thus providing an extension of the results 
presented by Reick 19|]. Then, we analyze the susceptibility describing the second harmonic 
response, providing a detailed discussion of K-K relations and of the sum rules, where 
the outputs of the simulations and the theoretical predictions are compared. Moreover, 
we define general, self-consistent algebraic relations connecting the harmonic generation 
susceptibilities of different observables and specialize them for the Lorenz 63 perturbed 
system treated here. As the Lorenz 63 system is far from being an Axiom A system, our 
results suggest that what shown in 2J] has a much wider applicability than what delimited 
by the hypotheses of the theoretical derivation. 

The paper is organized as follows. In Section 2 we recapitulate the general theoretical 
background behind our analysis. In Section 3 we describe our simulations, discuss our ex- 
perimental procedure, present the mathematical tools for data analysis, and obtain from the 
general theory specific results for the considered system. In Section 4 we describe our main 
findings. In Section 5 we summarize our results and present our conclusions. In Appendix 
lAl we present the formula for the first correction to the linear response. In Appendix [B] we 



show a procedure for defining self-consistent relation connecting the harmonic generation 
susceptibilities of various observables. 



II. THEORETICAL BACKGROUND 



We consider the dynamical system x = F{x). We then perturb the flow by adding a 



time-modulated component to the vector fie 
by X = F{x) + e(t)X{x). Following Ruelle 



d, so that the resulting dynamics is described 
, |2|], in the case of a perturbation to an Axiom 



A dynamical system, it possible to express the expectation value of a measurable observable 
$(a;) in terms of a perturbation series: 

/CXD 
Po(dx)<|. + 5^(<|.)(")(t) (1) 
n=l 

where po{dx) is the invariant measure of the unperturbed system. The n*'' term can be 
expressed as a n-uple convolution integral of the n*'* order Green function with n terms, 
each representing the suitably delayed time modulation of the perturbative vector field: 



oo cx) oo 

<^)^''\t)= J y ... y"daida2...da„Gi"^(ai,...,a„)e(t-ai)e(t-a2)...e(t-a„). (2) 



The n*'* order Green function G^^\ai, . . . ,an) is causal, i.e. its value is zero if any of the 
argument is non positive, and can be expressed as time dependent expectation value of an 
observable evaluated over the measure of the statistical mechanical system: 

G^^\ai, ...,an) = J po{dx) e(ai)e(cr2 - f^i) • • • ©(c^n - C^n-l) X 

xAn((T„ - an-i) . . . An(a2 - ai)An((Ti)<l>(x), (3) 

where 6 is the Heaviside function, A(«) = A(x)V(») describes the impacts of the pertur- 
bative vector field, and 11 is the time evolution operator due to unperturbed vector field so, 
that n(r)A(x) = A{x{t)) for any observable A. As discussed in 2J], in the case of quasi 
equilibrium Hamiltonian system po{dx) = po{x)dx corresponds to the usual (absolutely con- 
tinuous) canonical distribution, and the Green function given above is the same as that 



obtained in standard Kubo theory 



Instead, in the more general case analyzed by 



Ruelle, po{dx) is the singular SRB measure. 



If we compute the Fourier transform of we obtain: 



($)W| 



/ n \ 

J ■ ■ ■ J dt^i • • • du;nxi"^ (t^i, . . . , u;„) e(a;i) . . . e(cj„) x S lu -"^uji j 

-oo -oo ^ ' 



(4) 



where the Dirac 6 guarantees that the sum of the arguments of the Fourier transforms of 
the time modulation functions equals the argument of the Fourier transform of and 
the susceptibility function is defined as: 



dti . 



exp 



(5) 



Whereas the Green functions describe coupling processes in the time-domain, this function 
describes the impact of perturbations in the frequency- domain having frequency 'YTj=i^j- 

As discussed in l|, the linear susceptibility x^^^ (^) ^in analitic function in the up- 
per complex UJ plane, as a result of the causality of the linear Green function, so that 
it obeys Kramers-Kronig relations j^, lo|. 



22] . Following 



m 



241 ] it was shown 



that x^"* (tuu; cj, . . . , cij), responsible for n*"^ order harmonic generation processes, features 
the same analytic properties as the linear susceptibility. The asymptotic behavior of 
(ntu; UJ, . . . ,uj) is determined by the short-term response of the system, and it can be 
proved that in general [nuj; uj, . . . ,uj) ~ with uj oo, with /3 > and integer. 

It is possible also to prove that, if /3 + n is even, a = is real, whereas, if /3 -|- n is odd, 
a = iaj is imaginary. The values of a and (3 depend on the specific system under investiga- 
tion, on the considered observable, and on the order n of nonlinearity 24]- Therefore, since 
all moments of the susceptibility have the same analytic properties, a large set of pairs of 
independent dispersion relations is obtained: 







^^""Re{xf^ {nuj';uj',...,u;') 



oo /2p 

duj 



UJ'^ - UJ^ 



oo 

V^'Rejxf^ (na.;a;,...,u;)}=P J duj 



u"^+'lm{xi^ {nu;';u',...,u;')] 



UJ'^ - UJ^ 



(6) 



(7) 



with p = 0, . . . , 7 — 1, in order to ensure the convergence of the integrals, with P indicating 
integration in principal part. Comparing the asymptotic behavior with that obtained by 
applying the superconvergence theorem 2l|] to the general K-K relations (l6])-(I7j), we derive 



the following set of vanishing sum rules: 

oo 

y"cu'2PRe|xi"^ (W;o;',...,cu')}dcj' = 0, < p < 7 - 1, (8) 


00 

y"c^'2P+^Im|xf^W;c^',...,cu')}dcu' = 0, 0<p<7-2, (9) 


If /3 + n = 27 we obtain the following non-vanishing sum rule: 
00 

I u;"P+Hm{xtHnu'-u',...,u;')}du' = -an^, P = l-l, (10) 


whereas, if + = 27 — 1, the non- vanishing sum rule reads as follows: 
00 

Juj"^Re{x^^\nu;';u;',...,u;')}duj' = aj^, p = 7 - 1- (11) 


III. DATA AND METHODS 



Following Reick 19[ , we perturb the y equation of the classical Lorenz 63 system with a 



weak periodic forcing: 

X = a{y — x) 

y = rx — y — xz + 2e cos{ujt)x (12) 
z = xy — bz, 

where a = 10, 6 = 8/3, r = 28, u is the frequency of the forcing, and e is the parameter 
controlling its strength. Note that Reick used e instead of 2e; our choice will be motivated 
later. With these parameters values, the (unperturbed, e = 0) classical Lorenz 63 system 
is chaotic and nonhyperbolic, as changes in the value of r cause sequences of bifurcations, 
which alter the topology of the attractor. Nevertheless, we stick to the definition provided 
in Eq. [T]and define operationally the impact of the perturbation flow on the observable $(x) 
as follows: 

(5$(t, to, xo)) = <^ix{t, to, xo)) - ($)o (13) 

where xo and to are the initial conditions and the initial time, respectively, and (•)o = 
J po{dx)». We choose an initial condition belonging to attractor of the unperturbed system 



and, without loss of generality, we set to = 0. Since the system is chaotic, the Fourier 
Transform of 6^{t,xo) has a continuous background spectrum, so that, detecting and disen- 
tangling the response of the system to the cu-periodic perturbation, which results into peaks 
positioned at the harmonic frequencies mu, m > 1, is harder than in the quasi-equilibrium 
case. 



A. Linear and Nonlinear Susceptibility 

Since e{t) = 2ecos{ujt) = e {exp[iut] + exp[— icjt]), we obtain that the linear susceptibility 
can be expressed as: 

xS^(t^) = lim lim xi>i^^,xo,e,T) (14) 

where 

1 /-^ 

X$(a;,a;o,e,T) = - / dr(5$(r, xq)) exp[iu;r] (15) 

contains information on the full response (linear and nonlinear) of the system at the fre- 
quency of the forcing term. Note that, consistently with our definitions, there is a factor 



2 of difference with respect to the (equivalent) expression given in 19|. The susceptibility 
X^(^\itj), when limits are considered, does not depend on xq. Whereas the e — > limit corre- 
sponds to the physical condition of considering an infinitesimal perturbation, in the T ^ oo 
limit the noise-to-signal ratio goes to zero. Of course, in reality for every simulation we can 



measure x<s>{^,Xo,e,T), so that, as discussed in [l9|, T must be long enough for detecting 
any signal for finite time series. The larger the parameter e, the better is the signal-to-noise 
ratio, but, at the same time, the worse the validity of the linear approximation (and, in 
general, of the perturbative approach). 

As we wish to keep e as small as possible in actual simulations, and avoid using very long 
integrations, which may problematic in terms of the computer memory, in order to improve 
the signal detection we introduce an additional procedure based upon ergodic averaging. For 
each frequency component of the background continuous spectrum, which is related to the 
chaotic nature of the unperturbed dynamics, the phase of the signal depends (delicately) 
on the initial condition xq. Since, instead, the phase of the response of the system to 
the external poerturbation is (asymptotically) well-defined (see Eq. [T4l) . by averaging the 
susceptibility over an ensemble of initial conditions randomly chosen on the attractor of the 



unperturbed system we can improve the signal-to-noise ratio. We then define: 

1 ^ 

X^{u},e,T,K) = —^X'S>{^^,Xj,e,T) (16) 
i=i 

where xj are random initial conditions belonging to the attractor of the unperturbed system, 
such that: 

lim X^(,i^,(^,T,K) = / po(dx)x$(u;,x,e,T). (17) 

Therefore, we choose X'S>{^^^jT, K) as our best estimator of the true susceptibility x^^\^) 
in all of our calculations. 

We now propose some procedures for analyzing the nonlinear response of the system. 
Probably, the two most relevant nonlinear phenomena usually observed in weakly perturbed 
nonlinear systems are the correction to the linear response and the harmonic generation 
process. The correction to the linear response, which is basically described by the third 



order Kerr susceptibility, cannot be treated with the K-K formalism, as discussed in 10|, 



16 



so that it will not be discussed further. Anyway, an operational definition of such a 
susceptibility is provided in App. \M 

We then focus on the lowest order nonlinear susceptibility responsible for n'^^ order har- 
monic generation is x|^''(r2a;; cj, . . . ,uj) [l^, which, as previously described, obeys an entire 
class of K-K relations. Such a susceptibility gives the dominant contribution to the system 
response at frequency nuj, so that, along the lines of the formulas proposed for the linear 
susceptibility, we define: 

x'i\nu} ; u,..., u) =\im\im x$(na;, xq, e, T), (18) 



where 



1 

X$(ncj,a;o,e,T) = — / dr5<l>(r, Xq)) exp[iria;r]; (19) 
e Jo 



note that with Reick's parametrization of the periodic perturbation the factor before the 
integral would be 2"/e". Following the same approach described in Eqs. [T6lfT7l for the linear 
case, we choose x<s>{''^^j^jT, K) as our best estimator for the actual n*'^ order ra*'^ harmonic 
generation susceptibility. 

Reick 191] showed using symmetry arguments that if $(x, y, z) = x^y^z"^ with k + I odd, 
the linear response of the system to the periodic perturbation is vanishing. Therefore, as 
we wish to join on the analysis of linear and nonlinear susceptibility and improve the work 



by Rieck, we stick to his approach and concentrate on the observable ^{x,y,z) = z. In 
Appendix [Bl we show how results for harmonic generation susceptibility can be extended to 
general (polynomial) observables, thus extending the result obtained by Reick for the linear 
case. 

Concluding, we wish to specify that all integrations of the system shown in Eq. [T^] 
have been performed using a Runge-Kutta 4*^* order method with e = 0.25 up to a time 
T = 5000. For each of each of the considered values of u, which range from 0.0257r to 327r 
with step 0.057r, we have sampled the measure of the attractor of the unperturbed system 
by randomly choosing K = 100 initial conditions. We have analyzed the linear response at 
the same frequency u of the forcing plus the process of 2"'^ harmonic generation (n = 2). 



B. Asymptotic behavior, Kramers-Kronig relations, and sum rules 

The perturbation vector field introduced in the Lorenz 63 system [12] is X = (0, x, 0) and 
the modulating function is e(t) = 2ecos(u;t) = e(exp[iu;t] + exp[— icjt]). Following Eq. [3l 
the linear Green function G^^\t) results to be: 



Gi'H^) = J Po(dx)e(r)An(r)<l>(x(r)) = j po(dx)e(r)XV$(x(r)) = j po(dx)e(r)x9,^(r) 

(20) 

Since 

/cxD 
eit)t' exp[iujt] k\^^. (21) 
-CXD ^ 

the asymptotic behavior of the linear susceptibility is determined by the short term behavior 
of the linear Green function. We then Taylor-expand z{t) in Eq. [20] in powers of r by 
considering the unperturbed flow, integrate over po(dx), and seek the lowest-order non- 
vanishing term. We substitute z(r) = z + tz + o(r) in Eq. [3] and obtain: 

G^z\r) = j po{dx)e{T)xdyZ + T j po{dx)e{T)xdy{xy - hz) + o(r) = e(r)r {x")^ + o(r). 

(22) 

Since ^ 0, we have that: 



(23) 



which imphes that the real part dominates the asymptotic behavior, whereas the imaginary 
part decreases at least as fast as uj~^. This proves that the following K-K relations apply: 

--ImK'M}=p/d.'-l^ (24) 



^Re (xi^) (^)} = P / do;' i (25) 

2 J uj' — uj^ 







and that the following sum rules, obtained by considering the a; — >^ oo limit of the K-K 
relations, are obeyed: 

^i?i?« = 1 Re{x«(o;')}da;' = (26) 



oo 

SRI? = I c^'Im {x? [u')} du' = I (x^>„ , (27) 







where the pedix in the name refers to the moment considered in the integration. 

Regarding the second harmonic response, we have that by plugging the formula for the 
second order Green function presented in Eq. [3] in Eq. [5l setting uji = uj2 = u;, and substi- 
tuting Ti = 0"! and r2 = (72 — 0"!, we obtain the following expression for the susceptibility: 

oo oo 

xf (2cu; uj,uj)= j j dridr2e(ri)e(r2) exp [i2uti] exp [ic^r2] j po(dx)An(r2)An(ri)$(a;) = 

— oo — oo 
oo oo 

^ j j dTidr26(ri)6(r2) exp [i2uoti] exp [1UJT2] j po{dx)xdyX{Ti)dyz{Ti + T2). 
—00 —00 

(28) 

After a somewhat cumbersome calculation, performed along the lines of the linear case, we 
obtain that, asymptotically, 

X^'\2u;-u;,u;)r.a{x')Ju;\ (29) 

As in the linear case, the real part dominates the asymptotic behavior. Note that the 
calculation of asymptotic behavior for specific susceptibilities can be greatly eased by making 
use of the recurrence relations presented in Appendix [Bl 

We then obtain that the following set of generalized K-K relations are obeyed by the 



second harmonic susceptibility: 



(30) 



7 u''^^hm{xf\2u')\ 



(31) 



with p = 0, 1, where we have shghtly simphfied the notation. Note that in this case two 
independent pairs of dispersion relations can be established. When considering the a; — oo 
limit of the K-K relations, the following vanishing sum rules can be derived: 



oo 

SRR^o^ = J dco'Re {xi^^ (2cu')} = 



(32) 



= j dcu'cu'Im {xP (2cu')} = 



oo 

SRR^^ = J duo'uo'^Re {xf (2^.^')} = 0, 



whereas the only non-vanishing sum rule reads as follows: 



(33) 



(34) 



oo 

SRiP = I du'u^lm {xf^ {2u;')} = -^a {x\ 



(35) 



In actual applications in data analysis, both K-K relations and sum rules are negatively 
affected by the unavoidable finiteness of the spectral range considered. Whereas the zero- 
order solution to this problem is to truncate the integrals at a certain cutoff frequency 



^^cutoff, niore advanced processing techniques [l5|, l27|, l28| have been introduced in order to 
ease this problem, which, when only a rather limited frequency range can be explored, can 
greatly decrease the efficacy of the dispersion theory. 



IV. RESULTS 



A. Linear susceptibility 

In Fig. [T] we present the measured xi^'* {^) resulting from the definition and practical 
procedure described in Eqs. [T^lfTTl We observe as distinct features a peak in the imaginary 
part of the susceptibility and a corresponding dispersive structure in the real part. Good 
agreement with 19| is found. The resonant frequency 27ri^r with I'r ~ 0.9 roughly corresponds 
to the frequency dominating the autocorrelation spectrum of the z variable (not shown). 
The features of the susceptibility presented in Fig. [T] are in good qualitative agreement 
with the (proto)typical structures of the linear susceptibility of a forced damped Lorentz 
linear (or weakly nonlinear) oscillator model, having natural frequency ujq = 2'Ki'r fl^ - 
Nevertheless, note that in this case the amplification of the linear response does not result 
from a deterministic, mechanical resonance, and the autocorrelation spectrum does not 
give us an information wholly equivalent to that provided by the linear susceptibility. In 
fact, whereas the free fiuctuations of the system take place, by definition, on the unstable 
manifold, the forced fiuctuations, generically, do not obey such a constraint. This is the basic 
reason why the fiuctuation-dissipation relation does not apply in this context, as remarked 
in [241. 



The observed asymptotic behavior is in agreement with what obtained in Eq. [231 In Fig- 
[2] we show that for > 30 we have |xi"'^^(a;)| ~ {x^)q/uj'^ within a high degree of accuracy. 
This result allows for the actual application of the K-K relations given in Eqs. [2H Obviously, 
we deal with a finite frequency range, while Eqs. [2H require in principle an infinite domain 
of integration. Whereas more efficient integration schemes and data manipulation could be 



28| ]. or, more simply, the 



envisaged and implemented, such as those described in 
consideration of the asymptotic behavior, we simply truncate the dispersion integrals shown 
in Eqs. [24l at the maximum frequency considered in our simulations oJcutojf ~ 100. We 
present the result of such a simple Kramers-Kronig analysis for the linear susceptibility in 
Fig. [31 Note that, given the relatively slow asymptotic decrease, only one pair (p = 0, see 
pedix in the legend of Fig. [3]) of K-K gives convergence. We observe that the agreement 
between the reconstructed and the measured susceptibility is outstanding almost everywhere 
in the spectrum, with the only exception being the slight underestimation of the main 
spectral features, which is somewhat physiological as we are treating integral relations, 
which tend to smooth out the functions. 

In Fig. [Hwe plot the cumulative value of the integrals of the sum rules given in Eqs. 



up to Umax- The sum rule of the real part of the susceptibility converges to zero as predicted 
by the theory: extending the Umax range in Fig. H] and considering the asymptotic behavior, 
the curve approaches the x-axis. Moreover, the first moment of the imaginary part, which in 
the case of optical systems corresponds to the spectrally integrated absorption, converges to 
the predicted value ii/2{x^)q. This implies that, in principle, the measurement of the linear 
response of the system at all frequencies allows us to deduce the value of an observable of the 



unperturbed flow. This class o 



results are widely discussed and have been widely exploited 



in the optical literature 16|, l25|, where quasi-equilibrium systems are considered. 



B. Second harmonic generation susceptibility 

Rather positive results have been obtained also when considering the second harmonic 
generation susceptibility. The observed asymptotic behavior is in excellent agreement with 
what predicted in Eq. [29j in Fig. [2]we show that for > 20 we have |xi^'*(2c<j)| ~ (j{x'^)q/uj^. 
Therefore, in this case we try to verify the two pairs of independent K-K relations given in 
Eqs. [3OII3T] by performing the numerical integration along the lines described in the linear 
case. Results are presented in Fig. [6] for the real part and in Fig. [7] for the imaginary part 
of the susceptibility, where for simplicity only the main spectral features are presented. We 
observe that the agreement between the measured and K-K reconstructed susceptibility is 
quite good for both the p = and the p = 1 dispersion relations, except for the presence of 
overly smoothed spectral features and, in the p = 1 case, for the low frequency divergence. 
Such a pathology is typical of higher-order K-K dispersion relation, whose integrals feature 
a slower convergence, and can be cured by extending the investigated spectral range (see 



discussion in 



16 



231]). Anyway, these results are very encouraging and pave the way to 
the detailed analysis of the spectral properties of higher order response of perturbed chaotic 
systems. 

An old tale of nonlinear optics in solids says that, near resonance, the second harmonic 
susceptibility can in general be approximated according to the Miller's rule as follows: 

{2u) ~ Ai^^xi^^ (2u;) X^P {uf = X? {2u)mr , (36) 

where A]^J- is the so-called Miller's Delta, which is related to the nonlinearity of the potential 
of the specific crystal under investigation [l^. The Miller's rule, which can be extended also 
for higher-order harmonic generation processes, is related to the perturbative nature of 
harmonic generation processes [l^. We have constructed the Miller's rule approximation to 



the second harmonic susceptibihty, and resuhs are also shown in Figs. [MZl with AjJ ~ —1/3. 
The agreement is very surprising, and confirms that the functional form of the harmonic 
generation susceptibility is of very general character. Moreover, the Miller's rule can be very 
helpful in studying harmonic generation processes when only a limited amount of data on 
the nonlinear processes can be directly obtained from observations. 

The final step in this work has been the verification of the sum rules shown in Eqs. [3211351 
The first two sum rules are presented in Fig. [HI where we plot the cumulative value of the 
integrals in Eqs. [32l[33] up to ujmax- We observe that, as predicted by the theory, the integrals 
converge to with a rather satisfying accuracy as Umax goes to infinity. The other two sum 
rules given in Eqs. [34l[35| which are typical of the second harmonic generation susceptibility 
(no convergence is obtained in the linear case) are instead shown in Fig. [9] The second 
moment of the real part integrates to zero with high accuracy, whereas the third moment 
of the imaginary part converges to the value — 7r/2cr(a;^)o, as predicted by the theory. This 
confirms that a we have obtained a complete and self-consistent picture of second harmonic 
response. 

Note also that, by combining Eqs. [35] and [27], we can deduce the value of a, which is 
a fundamental parameter of the Lorenz 63 system, independently of the specific attractor 
properties (which are implicitly contained in the value of (a;^)o). This confirms that a 
detailed analysis of the response of the system to an external perturbation can in principle 
provide rather fundamental information on the structure of the unperturbed system. 



V. CONCLUSIONS 



This paper takes advantage of the scientific results obtained within various stream of re- 
search activities, such as the results by Ruelle on the linear and nonlinear response function 
for non equilibrium steady state chaotic systems [l, 2!, the related theoretical contributions 
of the author on the general theory of dispersion relations 2J], the numerical experimenta- 



tions by Reick [19|, the procedures for the spectral analysis of the linear [25] and nonlinear 
susceptibility in optical systems [igI, and provides the first complete analysis of the spec- 
tral properties of the linear and nonlinear response of a chaotic model - the celebrated, 
prototypical Lorenz 63 system - to an external, periodic perturbation. 

We have first provided a definition of the n^^ order (n > 1) harmonic generation suscep- 
tibility for the response of an generic observable to the external field perturbing the flow of 



the Lorenz 63 system and defined a practical procedure to extract it from the time series 



resulting from numerical integrations, which includes the computation of ensemble averages 
of the Fourier transform of the signal. Ensemble averaging drastically decreases the noise 
due to the internally generated variability of the system due to the chaotic behavior. A 
simple way for characterizing and computing the correction to the linear response due to 
finiteness of the perturbative forcing has been provided in App. [A] 

We have then shown how to theoretically compute ab initio, starting from the perturba- 
tive Green function, the asymptotic behavior of the linear and second harmonic generation 
susceptibility of the observable $ = 2;, and how to define the set of Kramers-Kronig relations 
and sum rules the susceptibilities have to obey. Moreover, in App. [HI have presented a set 
of general self-consistency relations allowing a for a extensive generalization of the obtained 
results. 

The results of the numerical simulations have shown that linear susceptibility features 
a strong resonance (peak for the imaginary part - dispersive structure for the real part) 
corresponding to the dominant frequency component of the z autocorrelation spectrum. The 
linear susceptibility, as predicted by the theory, decreases asymptotically as ~ —{x^)ouj~'^ 
and obeys K-K relations up to a very high degree of accuracy. Moreover, whereas the spectral 
integral of the real part of the susceptibility is found to be vanishing, the first moment of 
the imaginary part converges to 7r/2(x^)o, in agreement with the theoretical results. These 
results extend and give a more solid framework to the findings by Reick fiol . 

The results of the numerical simulations are in excellent agreement with the theory de- 
veloped also when the second harmonic process is considered. The second harmonic suscep- 
tibility decreases asymptotically ~ cr(x^)oC(;~^, and obeys precisely K-K relations. Moreover, 
thanks to the fast asymptotic decrease, K-K relations can be established for and are obeyed 
by the second moment of the second harmonic susceptibility. The resulting sum rules - ze- 
roth and second moment of the real part, first and third moment of the imaginary part are 
precisely obeyed by the measured susceptibility, with the former three being vanishing, and 
the last one converging to — 7r/2cr(x^)o. Therefore, by combining the two non vanishing sum 
rules for linear and second harmonic susceptibility, one could in principle not only obtain 
information on the expectation value of an observable (specifically, x^) at the unperturbed 
state, but also the on value of one of the parameters of the Lorenz 63 system (a). 

This is the first detailed theoretical analysis and numerical evidence of the validity of 
Kramers-Kronig relations and sum rules for the linear and the second harmonic response 
for the perturbed Lorenz 63 system. Our findings confirm the conceptual validity and 
applicability of the linear nonlinear response theory developed by Ruelle. As the rigorous 



theory basically deals with perturbations to an Axiom A system featuring an attractor with 
an SRB measure, the procedures defined and the results obtained here have to ben considered 
as heuristic, but, possibly, robust enough to be indicative of future research activities, both 
theoretically and numerically oriented. 

Since the fiuctuation-dissipation relation cannot be used straightforwardly in non equi- 



librium steady state systems, in spite o, 
examples in fiuidodynamics and climate 



several attempts in this direction {e.g. see some 



30|), due to non equivalence between forced and 



free fiuctuations, dispersion relations, which are based only upon the principle of causality, 
should be seriously considered as tools for diagnosing the behavior of these systems, under- 
standing their statistical properties, and reconstructing their output signals. In particular, 
dispersion relations could play a role for interpreting climate change signals and understand- 
ing in greater depth the dynamical processes involved in climate change 3l| . As a step in this 
direction, the author foresees adopting the mathematical tools presented in this paper for 
analyzing the response of a simplified climate model 32] to various kinds of perturbations. 
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APPENDIX A: CORRECTION TO THE LINEAR RESPONSE: KERR TERM 



The lowest order correction to the hnear response of the system at frequency u; can be 



expressed as a third-order term, responsible for what in optics is called the Kerr effect [16 1 
The corresponding nonlinear susceptibility x^i\^^ —lu,u) can operationally be defined as: 



x'i\uj;u;,-uj,uj) = \im\im ^ a;o, e, T) - x$ 



UJ 



(Al) 



Since the e — > limit is not attainable in reality, a good approximation to the Kerr suscep- 
tibility can be expressed as follows: 



x'i\uj, -Lj,uj) ~ \ {x<s>{uj,ei,T,K) - x<i>{uJ,e2,T, K)) , 
^1 



(A2) 



with 62 ^ ci. This term is responsible for the difficult y en countered by Reick 
an e-independent linear susceptibility. As well known 



10 



15 



16 



n 



191 ] in defining 



18| , the Kerr term does not 



obey K-K relations, and other techniques, such as maximum entropy method, are needed 
to obtain the real part of the susceptibility from the imaginary part or viceversa [l^. This 
has been confirmed in the present investigation (not shown), where the Kerr susceptibility 
has been constructed via Eq. IA2I using ei = 1 and €2 = 0.25. 



APPENDIX B: RESPONSE FUNCTION FOR GENERAL OBSERVABLES 



If we perform integration by parts in eqs. [T8lfT9| we obtain that: 



niujx'i\nuj;uj,. . .,uj) = x^^\nu;;uj, . . . ,u;), (Bl) 



where the hnear case is obtained by setting n = 1. Let's consider $ = x^y^z^ and substitute 
the hmit presented in Eqs. [T8]fT9l for the second member: 

1 d 

—niuoXklmi^^'^^^---^^) ~ ^™ / dr— (5[a;^y'z'^(r, xo)] exp[inci;r] = 

1 

= hm hm — / dT{k[F^ + e(t)Xj5 \x^-^y^ z'^ir, xq)] + l[Fy + e{t)Xy\5 \x^y^-^ z'^ij, xo 
+ m[F^ + e{t)X;,]S [x^?/'2"'"^(r, Xq)] } exp[inu;r] 

1 

= lim hrn— / dr {A;cr(5[x''"^y'+^2™(r, xq)] - kaS[x''y^z"'{T,Xo)]} exp[mujT] + 

1 

+ \im lim — y dr {/r5[x^'+^y'~^z™(r, Xq)] - W[x''?/'z"'(r, Xq)] } exp[irau;r] + 

+ hmhm^ / dr{-/(5[x'=+^?/'"^z™+^(r,Xo)] +m(5[x^+^2/'+^z™-^(r,Xo)]}exp[incjr] + 

1 /"^ 

+ hm hm— / dT\^—mbS[x^y''z"^{T,Xo)]^exp[inujT] + 

T oo e 6 

+ hm hm— / dr {/e (exp[ici;t] + exp[-icut]) ^fx^'+^y' ^2™(r, Xq)]} exp[inci;r]. 



(B2) 



where the pedices k,l,m refer to the monomial x^y^z"^ defining the observable. In Eq. 
we have used explicitly the perturbed Lorenz 63 flow given in Eq. [12] in order to compute 
the derivatives. The last term of the summation in the second member of Eq. IB2I can be 
rewritten as: 



1 



lim lim / dr {/5[x''+^y'"^^"'(r, Xq)]} exp[i(n + l)cjr] + 

, r,.r .^1 ,-^ rr., . 



hm \im-^ dr{W[x'=+^?/'-^^"^(r,xo)]}exp[i(n- 1) 

The first term is vanishing (the value of the integral ~ e'"+^ as it corresponds to an n + 1*'' 
order harmonic generation process), whereas the second term is just / times the susceptibility 
descriptive of the process of — 1*^ harmonic generation for the observable $ = x'^~^^y''~^z"^. 



Substituting the definition given in Eq. [3] in Eq. IB21 we obtain: 



- lxtii,i-i,m+ii'^) + mxtii^i+i^^.iinu;) + lxtti,i-i,miin - l)w), 

(B4) 



191 ] at all orders of nonlinearity. Obvi- 



which generalizes the n = 1 case presented in Reick 
ously, this procedure can be generalized to any dynamical system perturbed with any sort 
of periodic perturbation field, and allows for defining several self-consistency relations (see 
the second line of Eq. IB2p . 

As Eq. IBll can be easily generalized to all orders of derivatives k as 
{—niu)"^ X^^\nuj; u, . . . ,uj) = X^"/dtfe{$}('^'^' ■ ■ ■ possible to define a class of con- 

sistency relations analogous to that presented in Eq. IB4[ Nevertheless, one must be aware 
that the number of terms involved in these relations grows exponentially with the degree k 
of the derivative considered. 




FIG. 2: Asymptotic behavior of the absolute value of the linear susceptibility. The agreement with 
the theoretical prediction ~ {x^)qU}~'^ is apparent. 




FIG. 3: Reconstruction of the linear susceptibility with Kramers-Kronig relations. The agree- 
ment with the measured susceptibility presented in Fig. [T] is remarkable, except for the slight 
underestimation of the localized spectral features. 




FIG. 4: Cumulative value of the sum rules given in Eqs. [2611271 with upper integration limit set 
to 

^max- The extrapolation to infinity is in excellent agreement with the theoretically predicted 
values. 




FIG. 5: Asymptotic behavior of the absolute value of the second harmonic susceptibility. The 
agreement with the theoretical prediction ~ a{x'^)oUJ~^ is apparent. 
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FIG. 6: Real part of the second harmonic generation susceptibility in the resonance region. From 
left to right: measured susceptibility, reconstructed susceptibility via K-K relations with a = 0, 
reconstructed susceptibility via K-K relations with a = 1, reconstructed susceptibility via Miller's 
Rule approach. Agreement is remarkable 
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FIG. 7: Imaginary part of the second harmonic generation susceptibihty in the resonance region. 
From left to right: measured susceptibihty, reconstructed susceptibihty via K-K relations with 
a = 0, reconstructed susceptibility via K-K relations with a = 1, reconstructed susceptibility via 
Miller's Rule approach. Agreement is remarkable 
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FIG. 8: Cumulative value of the first two sum rules given in Eqs. [32]|33] with upper integration 
limit set to oomax- The extrapolation to infinity is in excellent agreement with the predicted values. 



400 




FIG. 9: Cumulative value of the last vanishing sum rules given in Eq. [34] and of the non- vanishing 
sum rule given in Eq. [35] with upper integration limit set to tdmax- The extrapolation to infinity is 
in excellent agreement with the predicted values. 
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